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O . Abstract: Non-perturbative production of fermions after chaotic inflation has been C_]) 

^ ! the object of several studies in the very recent past. However, the results in the most 
interesting case of production of massive fermions in an expanding Universe were so 

Ph! far known only numerically. We provide very simple and readable analytical formulae, v.^ ^ 

Oh! both for the spectra of the created fermions and for their total energy density. Their 

r-| ! derivation is closely related to the one adopted for bosons and exploits the fact 



that the production occurs during very short intervals of non-adiabatical change of 
^ ! the fermionic frequency. Our formulae show the presence of resonance bands if the 

expansion of the Universe is neglected, and their disappearance when the latter is 
included. As in the bosonic case, this last effect is due to the stochastic character 
that the expansion gives to the whole process. Backreaction is considered in the 
last part of the work. All our analytical results are in excellent agreement with the 
previous numerical ones in the regime of validity of the latter. However, a more 
accurate scaling for the energy density of the produced fermions is here found. 
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1. Introduction 

One of the key ingredients for our understanding of the early Universe is the mecha- 
nism of inflation |[l|, which constitutes a very elegant solution to several cosmological 
problems. 

Despite of the simplicity of the general idea, the details of the physical processes 
which govern it are still somehow unclear and matter of intense work. The two main 
aims of these studies are 

(i) to embed inflation in a context more motivated by particle physics (for a review 
see 0) and 

(ii) to understand the reheating phase which converts the energy density that 
drives inflation into the matter and radiation that we presently see. 

This second issue has been deeply influenced in the last decade by the possi- 
bility of particle creation through parametric resonance . The application of this 
phenomenon to creation of matter after inflation has been called preheating in the 
paper g], since (with the exception of some very recent versions 0) it is usually 
followed by a stage of (ordinary) perturbative reheating. 



Preheating of bosons 0, is characterized by a very efficient and explosive 
creation, due to the coherent effect of the oscillations of the inflaton field. This allows 
significant production even when single particle decay is kinematically forbidden. It 
has been very recently noticed that preheating of fermions can also be very efficient 
despite the production is in this case limited by Pauli blocking. 

Parametric creation of spin 1/2 fermions has been the subject of some works 
in the past. Pure gravitational production has been examined in refs. while 



creation by an oscillating background field is instead considered in the works |]10 



111 , |T2|, ^ 0, |T^. References |jTO|, |TT[ report results for creation in a Minkowski 
space. Reference |]12| studies the production of massless fermions after a X(j)^ inflation, 
exploiting the fact that this case can also be reconducted to a static one. In this work, 
production in a static Universe after chaotic inflation is also considered, and some 
conjectures on the effects of the expansion are made. Moreover, the full calculation 
of preheating of massive fermions after chaotic inflation in an expanding Universe 
has been performed numerically in ref. |^. 

These last works had a great impact on the most recent studies. For exam- 
ple, their results turned also useful to the study of gravitinos production at pre- 
heating [14|. This issue is particularly important, since gravitinos can easily over- 
close the Universe (if they are stable) or (if they decay) spoil the successful pre- 
dictions of primordial nucleosynthesis through photodissociation of the light ele- 
ments. Gravitinos can be thermally produced during the stage of reheating. To 
avoid this overproduction, the reheating temperature Trh after inflation cannot be 
larger than ~ (10^ — 10^) GeV |jl5[. However, it has been realized [0] that the 
non-thermal production of helicity ±1/2 gravitinos (whose equation of motion can 
be reconducted to the one of an ordinary spin 1/2 Dirac particle) can easily be 
more efficient than the thermal one, and this in general leads to more stringent 
upper bounds on Trh- Several papers related to the works [l^ have recently ap- 
peared [|l^, |l^ . 

Another important implication of preheating of fermions is constituted by lepto- 
genesis, as the work and the related papers |]18| show. In this scheme [jl9[, a lep- 
tonic asymmetry is ffist created from the decay of right-handed neutrinos, and then 
partially converted to baryon asymmetry through sphaleronic interactions. Since lep- 
togenesis is very sensitive to both the mechanism of creation of the heavy neutrinos 
and to the neutrino mass matrices, it could constitute an interesting link between 
preheating and the experiments on neutrino oscillations. 

Other phenomenological implications of these works appear in refs. [|^, with 
preheating as a possible mechanism for creating superheavy relic particles responsible 
for the ultrahigh energy cosmic rays, and in ref. |21], where the possible impact of 
fermions produced during inflation on the microwave background anisotropics and 
on the large-structure surveys is considered. Finally, fermionic production can play 
an interesting role in hybrid inflationary models |T^, 



Due to the large number of these studies, it may be worth to reconsider the 
basic mechanism of fermionic preheating. In particular, it should be important to 
give an analytical confirmation to the results of production of massive fermions in the 
expanding Universe, which are so far known only numerically from the analysis [0]. 
This is the aim of the present work. 

In the next section we revise the basic formalism for preheating of fermions. We 
consider creation of very massive particles right after chaotic inflation. The coupling 
of the fermions to the oscillating infiaton gives them a time varying mass. As it is 
known, this can cause a non-adiabatical change of the frequency of the fermions and 
their consequent creation. In case of very massive fermions, the non-adiabaticity 
condition can be satisfied only when their total mass vanishes, and the production 
occurs at discrete intervals, until the infiaton oscillations become too small for the 
total fermionic mass to vanish. 

In section ^ we derive analytical formulae for the spectra of the fermions after a 
generic production. Our derivation follows the one developed in ref. [^] for preheat- 
ing of bosons. It exploits the fact that the production occurs in very short intervals 
around the zeros of the total fermionic mass: the calculation is made possible from 
the fact that the occupation number can be considered as constant outside these 
small regions, and that the expansion of the Universe can be neglected inside them. 
As a result, the only physical quantities relevant for the creation are the time deriva- 
tive 0' of the infiaton field and the value of the scale factor a at each production. 
The derivation of the "fermionic counterpart" of the formulae obtained for preheat- 
ing of bosons in ref. |^ has also been done to a certain extent in the work [|21]], where 
the results of a single production during inflation is given. However, when one is 
interested in the successive productions, a more detailed study is necessary, as our 
analysis shows. As may be expected, the final results that we obtain closely resemble 
the ones of the work What is most surprising is their excellent agreement with 
the numerical results, as some figures provided manifest. 

In section H we consider the production in a non-expanding Universe. In this case 



our analytical formulae considerably simplify and agree with the ones of refs. [11, 12 



In particular, they show the presence of resonance bands which are anyhow limited 
by Pauli blocking. 

In section ^ we study the more interesting case of production in an expanding 
Universe. As indicated in ref. and as confirmed by the numerical results of 
the work [0], the creation is now very different with respect to the previous case. 
The expansion removes the resonance bands and the production (almost) saturates 
a Fermi sphere up to a maximal momentum. Our analytical results confirm this 
behavior. 

In section ^ we also calculate the total energy density px of produced fermions, 
which may be the quantity of most physical relevance. We compare our results with 
the ones of ref. where it is shown that the final value of px (normalized to the 



inflaton energy density) scales linearly with the parameter q = ((7^0Q)/(4m^)/ while 
it depends very weakly on the fermion bare mass mx- However, these numerical 
results are valid only in a limited range for mx, and it has been wondered if the 
density px decreases at values of mx below this range. Our analytical results give 
a positive answer in this regard. To see this, a proper average of the analytical 
formulae must be done, exploiting the fact that the expansion of the Universe gives 
the production a stochastic character. In this way one can get a "mean" function 
that interpolates very well between the maxima and the minima of the spectra of 
produced particles. Again we derived it in close analogy with what is done in the 
bosonic case, and again the results that we get are in very good agreement with the 
numerical ones of ref. in the region of validity of the latter. 

All this analysis neglects the backreaction of the produced fermions on the evo- 
lution of the inflaton field and of the scale factor. Despite the difficulty of a more 
complete treatment, backreaction effects can be understood at least in the Hartree 
approximation. This was done numerically in ref. 0]. In section |^ we see that the 
analytical formulae here provided allow to understand the effects of backreaction 
observed in the numerical simulations. 

2. Production of fermions at preheating 

In this section we revise the basic formalism for production of fermions at preheating. 
Our presentation follows the ones of refs. 0, |2l|] with the correction of few typos. 

We start from the Dirac equation (in conformal time 77) for a fermionic field X 
in the Friedmann-Robertson- Walker background: 

(^Wd^, + t^H^''-m^X = 0, (2.1) 

where a is the scale factor of the Universe, H = a'/a"^ the Hubble rate,^ and m the 
total mass of the fermion. 

Fermionic production is possible if this last quantity varies non adiabatically with 
time. This may happen during the coherent oscillations of the inflaton (p at reheating 
in presence of a Yukawa interaction (pXX, such that the total mass is given by 

m{r]) = mx + g(l){r]) ■ (2.2) 

We will always consider fermions with very high bare mass mx- In this case, non 
adiabaticity can be achieved if the Yukawa interaction is sufficiently strong to make 
the total mass ( p.2|) vanish. Fermions are then created whenever the inflaton fleld 
crosses the value 0* = —mx/g- 

^In this expression g is the Yukawa couphng between the inflaton and the fermions, 0o the initial 
value of the inflaton, and its mass. 

^Here and in the following prime denotes derivative which respect to ?/. 



We will study the production after chaotic inflation, that is while the inflaton 
field coherently oscillates about the minimum of the potential 

V=^ml<p\ m^-lQ^^GeY. (2.3) 

If one neglects the backreaction of the created particles (this effect will be considered 
in the last part of the work), then, after few oscillations, the inflaton evolves according 
to the formula 

0(t)^^^^^, (2.4) 



where t is the physical time. The presence of the t at the denominator of eq. ( |2.4| ) 
shows the damp of the oscillations due to the expansion of the Universe. It thus 
follows that there exists a final time after which |0| < Tnx/g and the total mass no 
longer vanishes, so that the production ends. 

To proceed with our analysis, we redefine x = a'^^'^X. Equation (p.l|) acquires 
the more familiar form 

- am) X = . (2.5) 

We decompose 

^^-^e-^- J2 [Mk, V)ar{k) + vr{k, vHi-k)] , (2.6) 

with^ 

Vr{k) = CuJ i-k) . (2.7) 

It is convenient to take the momentum along the third direction k = kz and 
to define 



Ur = 



T 



Vr = 



T 



:2.8) 



with ijjj^. = ^^j and ■?/'_ = (^^j eigenvectors of the helicity operator cr ■ v/|v|. From 
the normalization adopted in eq. (|2.8| ), one can impose 



u^Us = v^Vs = 5rs ■, u^Vs = Q-, |m+P + |m_P = 2. (2.9) 
With this choice, the Dirac equation ( p.5|) rewrites 

u'j^ijf) = iku^irf) =1= iamu±{ri) , (2.10) 



^We choose 



1 0\ 1 _ ( ^ -i(T2\ 2 _ ( ^ ^f^l 



^ ' -ly ' ^ \-i(T2 y ' ^ 
^ l-i Ui 



which can be decoupled into^ 

ul + [iul ± i{amy] u± = , uj^{r]) = k"^ + a'^m'^ . (2.11) 
Under this decomposition, the hamiltonian of the system 

H=lj d'^x^ix) i-id,) x{x) (2.12) 

rewrites 

H=^ f d^kJ2 {Ek iv) K (^) (^) - (k) hi {k)] + 

+ (r/) K i-k) ar (k) + F; (t]) aj (k) bl (- A;) } , (2.13) 
where (using the equations of motion) 

Ek = kRe + a m (l — u'^u+) , 

k 

Fk = -{ul-ul) +amu+u-, ^fc + l^fcP = ^^- (2-14) 

It is always possible to choose an initial configuration with the hamiltonian in 
the standard (diagonal) form. If we take, for example, 

1 =F j G , arbitrary phase , (2-15) 

we have Ek{r]o) = u, Ffc(?7o) = 0. 

However, the evolution equations ( p.ll|) drive F^ different from zero and a di- 
agonal form for H can be recovered only after a Bogolyubov transformation on the 
creation/annihilation operators 

a(?7, k) = ak{r])a{k) + (3k{y])h'^ {-k) , 
6+(r/, k) ^ -mv)a{k) + al{r^)b+{-k) . (2.16) 

While it is immediate to see that the equal time anticommutation relations on 
both sets {a^'^\b^^^} and {d^'^\b^'^^} enforce 

|afc|' + |/3fc|' = l, (2.17) 

the use of some algebra shows that a diagonalization of the hamiltonian is possible 
with the choice 

C-P(^), l/^P^^J^^^. (2.18) 



F* y' 2uj{E + iu) 2oj 



^We do not deal with the analogous equations for v±^ since from condition (^j) it simply follows 

v+{ti) = -u*_{T]),v_{r]) = u*^_{r|). 
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The time varying operators a^^-* and 6*^+-* are employed to define time dependent 
Fock spaces, each of them built from the zero (quasi)particle state at the time rj: 

a(77)|0,) = %)|0,)=0. (2.19) 

Let us consider the state |0^p) = |0) which describes a system with zero particles 
at the initial time rjo. This is no longer true for generic time rj > rjo, since the 
occupation number at t] is defined in terms of the operators which diagonalize the 
hamiltonian at that moment. The particle density per physical volume \^ = at 
time 7] is indeed given by^ 

niv) ^ (0|^|0) = {0\^Y1 I 7^a.^(^,r)a.(r/,r)|0) = 

= ^ [ dkP\Pk\\ (2.20) 

which is different from zero whenever the hamiltonian is non diagonal in terms of 
the initial operators. The occupation number of created fermions is thus given by 
nk = |/5fep and the condition (|2.17| ) ensures that the Pauli limit < 1 is respected. 



3. Analytical evaluation of the occupation number 

In this section we calculate analytically the evolution of the Bogolyubov coeffi- 
cients ( |2.16D during the oscillations of the inflaton field after chaotic inflation. 



In this derivation, we exploit the fact that, in the regime of very massive fermions 
we are interested in, the creation occurs only for very short intervals about the 
points 0* = —mx/g where the total fermionic mass (cfr. eq. (|2.2|) ) vanishes. As the 
perfect agreement with the numerical results will confirm, this consideration allows 
one to treat the fermionic production with the same formalism adopted in the bosonic 
case [0. While far from the zeros of the total mass m the Bogolyubov coefficients 
can be treated as constant, whenever crosses (p^, a sudden variation occurs. Since 
the interval of production is very narrow, one can safely neglect the expansion of 
the Universe during the production and also linearize the function ~ 0^, + 

(p'{v*){v ~ V*)- As a consequence, the frequency lu defined in eq. ( |2.11|) acquires 
the form 

uj^^e + [r],) [r],) (v - r],)\ (3. 1) 

and the whole calculation strongly resembles the one for scattering through a qua- 
dratic potential. 

The use of this formalism is very well established in case of production of 
bosons @]. For what concerns fermions, it has been recently adopted by Chung 



et al. for the study of the production during inflation. Fermionic production 



'Of course the same result is achieved also for antifermions. 



during inflation is possible only if the coefficient g of the Yukawa interaction has op- 
posite sign with respect to the inflaton field during inflation, or the total mass (|2.2| ) 
would never vanish.^ If this is the case, it is possible to choose the value of g such 
that the production occurs only once during inflation, while during reheating |0| is 
always too small for having creation. Our derivation is strongly inspired by this 
work. However, we are interested in couplings for which the production occurs sev- 
eral times during reheating, and it is not at all guaranteed a priori that an analytical 
approximation may work also in this case. The present analysis not only positively 
answers to this question, but also provides very simple and explicit formulae valid 
for arbitrary number of productions. 

As may be expected, the final formulae are very similar to the ones found for 
production of bosons 0]. What is most surprising is the almost perfect agreement 
between the analytical results and the numerical simulations, as we will show with 
some figures at the end of this section. 

The first step for the derivation of the analytical formulae is to consider asymp- 
totic solutions of eqs. ( p.lO[ ) and ( |2.11| ). We look for solutions valid for (p not very 



close to 0*, where the adiabaticity condition u' <^ uj"^ holds. 

In this regime, the most general solutions of eqs. (|2.10 ) and are 



b[i- ^-^Y' e-r^.^. + a{i+ "-^y" e^r^.^. , (3.2) 



with |y4p + \B\^ = 1 following from the condition + = 2. 

We can always put solutions of eqs. ( p. 101 ) and ( |2.11| ) into the form ( |3.2| ), with 



A and B in general functions of 1]. However, in most of the evolution (whenever the 
adiabaticity condition holds) it is a very good approximation to treat the coefficients 
A and B as constant. 

Substituting the expressions (|3^ ) into eqs. (|2.14|) , one finds 



F = 2ujAB, E = - , (3.3) 

from which it follows 

One can thus choose (up to an irrelevant global phase) 

A = a — ^B = p*, (3.5) 

where a and f3 are nothing but the Bogolyubov coefficients we are interested in. 

^Of course this constraint does not apply to our case, since the inflaton field changes sign after 
each half oscillation. 
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Notice that the initial choice A{rio) = 1, B{r]o) = corresponds to the initial 
condition (|2.15|) and to the zero particles state chosen in the previous section. 

As we have said, these coefficients undergo a sudden change whenever crosses 0o 
and then they stabilize to new (almost) constant values. Our aim is to find the values 
at the end of the variation in terms of the ones prior to it. 

We have not specified the lower limit of the integrals appearing in eqs. ( |3.2| ). For 
present convenience we choose it to be the time 77*1 of the first production (that is 
when = 0* for the first time). 

Let us consider the evolution equation ( |2.11| ) near the point 77*1. Since for high 
mass rrix the fermionic production is limited to a very short interval, one can neglect 
the expansion of the Universe during it and write the equation for 0(77) in a linearized 
form. We can thus write 

am(r7) ~ a*i5(0'^J?7 - r^.i) , a^i=a{r]^i), 0,i = ^ . (3.6) 

dr] v*i 



Following the notation of |21[], we define 

k 



P 



A/fi-l^^ila*! 



In terms of these new quantities, eqs. (|2.11|) rewrite (dot denotes derivative with 
respect to tY 

u± + {p^Ti + r^) M± = . (3.9) 

The point 17*1 is thus mapped into the origin of r and the region of asymptotic 
adiabaticity is at large |r|. 

In the asymptotic solutions ( |3.2|) we can see the behaviors 



am\i/2 p 
1 + 



00 J y/2T 

1 ^V2, 



UJ / 



±ipV2 



^For the production at the moment 77,„, when the total mass vanishes for the n-th time, eq. (3.9) 
must be replaced by 

u± + (/ ± z sign (0'„J + r") w± = . (3.8) 
The effect of this replacement on the final results is reported below. 
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for r — >• +00, and 



1 + 
1 - 



am \ 1/2 



am\ 1/2 

LjJ ) 



V2, 



P 



T 



P 

-2t 



(3.11) 



for r ^ — cxo. 

Equations (p.9|) are solved by parabolic cylinder functions D\{z) |23|. More 
precisely, the combination that matches the asymptotic solution ( |3.2|) at r — > — 00 is® 



u+{r) = A 



V2( 



p_ 
\V2 



l+jp2/2 



e^(f-^)e-V/8^_^_ 



+ r) + 



« r 



(3.12) 



In the above expression, A~ and 5^ denote the values of the Bogolyubov co- 
efficients before crosses 0^,, while the functions which multiply them are exact 
solutions of the linearized equation (p.9|) . The analytical approximation consists in 
considering them as solutions of the true evolution equation ( |2.11|) . For r +00 it is 
convenient to rewrite the solution ( |3.12 ) in terms of two different parabolic cylinder 
functions^ 



r(i + ip2/2) 



X 



X 



+ 



X 



V2( 



p_ 
\V2 



-ip2/2 



-A-e 



- np^ /A 



r(l-2p2/2) 



V2( 



^ ^'^^^ ^'e^(^-^)e--^^/«D_i_,,./2((l +z)r) 



(3.13) 



In this new expression, the functions within curly brackets correspond to the 
asymptotic forms at r = +00 of the two terms of the solution ( p.2| ). The coefficients 
in front of them give thus the new Bogolyubov coefficients in terms of the old ones. 

All this derivation can be easily generalized when productions at successive zeros 
of the total mass m are considered. The only important points are 

^We deal only with the function u+, since the study of u_ leads to the same results. 

^This rewriting is always possible since eq. ( |3.9D has only two linearly independent solutions. 
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(i) a difference in the values of the scale factor a and of the derivative (f)' at differ- 
ent T^^i's, 

(ii) a change of sign in the transfer matrix (the one which gives the new coefficients 
in terms of the old ones) whenever cj) crosses (p^, from below to above (cfr. the 
footnote just before eq. (|3.9|)), and 

(iii) the phase e^*^'^*'^'' which accumulates between rj^i and the rj^i considered. 
Putting all this together, one has 



for n odd . 



I^nj V -H* F* J \ I3l_^ 
Hn < — > —Hn for n even , (3-14) 



where an, Pn are the values of the Bogolyubov coefficients after the n-th production, 
and where 

Fn = V^l - e-7rp2gi(-f-argr(ip2/2)+4lii(Pn/2)-p^/2) ^ 



We remind that p„ = kj a/^I^Qo^. 

If one starts with no fermions at the beginning, one may choose ao = Ij /^o = 0. 
Then, applying successive times the "transfer" matrix ( |3.14 ), one can get the spec- 
trum of fermions produced after every 

Of course our calculation reproduces the result 

iVi = = e-"Pi (3.16) 

reported in ref. |^ . 



We numerically integrated the evolution equations for u±, 0, and the scale 
factor a}^ The results obtained with the analytical expression ( |3.14|) are always 
in very good agreement with the numerical ones. Just to give a couple of examples, 
we present here two cases at different regimes (we show them only with illustrative 
purpose, and the values of the parameters chosen have no particular importance). 
In figure |1| we present the spectrum of the fermions after two productions, that 
is after one complete oscillation of the inflaton field. In analogy with the bosonic 
case, we measure the strength of the coupling inflaton-fermions with the quantity 
q = 5f^0Q/(4m^), where 0o — 0.28Mp is the value of the inflaton at the beginning of 
reheating. In figure we choose q = 10^, while we fix the bare fermion mass to be 
mx = lOOm^. In figure ^ we show instead the resulting spectrum after 7 productions 
in the case q = 10^, mx = 4:171^. 



"'^'^Our starting point is at (/)(0) = 0.28A/p, short after inflation, <j>'{Q) = — O.lSMpm^ (as follows 
from a numerical evaluation of the inflaton alone during inflation), and a(0) = 1. 



11 



1 



0.01 



0.001 



0.0001 



1e-05 



1e-06 



1e-07 



analytical 
numerical 



k/ hih. 



100 



Figure 1: Spectrum of the fermions after two productions for q = 10^ and nix = 100m 
The expansion of the Universe is taken into account. 
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Figure 2: Spectrum of the fermions after seven productions for q = 10^ and mx = 
The expansion of the Universe is taken into account. 
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4. Production in a non-expanding Universe 

In the bosonic case, the study of the non-perturbative production in a non-expanding 
Universe has proven very useful in understanding the effects of the production. It 
is shown in ref. that the bosonic wave function satisfies the Mathieu equations, 
whose solutions are characterized by resonance bands (in momentum space) of very 
"explosive" and efficient production. It is then shown that, due to the expansion of 
the Universe, modes of a given comoving momentum k cross several resonance bands 
during the evolution. This gives the creation the stochastic character described in 



the work In ref. |T2| it is understood that an analogous behavior occurs also for 
fermions. The expansion is expected also in this case to spoil the clear picture of 
distinct resonance bands. This fact may help the transfer of energy to fermions, since 
the resonance bands in the fermionic case are anyhow limited by the Pauli principle. 
The expansion allows thus new modes to be occupied, and the production is no longer 
limited to the regions of resonance. In ref. fl^ it is stated that the production should 



''max- 



then almost completely fill the whole Fermi sphere up to a maximal momentum 
This behavior is confirmed by the numerical results of the work . In this section 
we will see that our analytical formulae can reproduce the resonance bands, while in 
the next one we will discuss the effects of the expansion of the Universe. 

Let us consider the matrices (we drop the index n in the matrix elements since 
all the pnS have now the same value p) 

F G\ /e^'^? \ /e^'^i 

-G F*J ' \ e-<) ' ^~ V e-'^ 

with G = e''^'p'^l'^ and 'd\ = p*^ uj^dr] {F was defined in the previous section). 

Without the expansion of the Universe, the inflaton field has the periodic 
evolution 

(j){ri = t) = (po cos {m^T]), (4.2) 
and all the -t^l ( |3.14| ) are hence sums of -^l and 'd^ (remember the rj^i are the moments 



at which the total fermionic mass vanishes). 

After the generic n-th complete oscillation one thus has 



q2n + l 



P;j V ^^.^^^^ jT2M'T,M...T,M'T,Mi^^j (4.3) 

with the combination O = T2M'^TiM repeated n times. 

One has thus simply to study the eigenvalue problem for O (notice detO = 1). 
This operator has the form 

o=(-* ^ 

V -B* A* 

B = FGei''"^''") - F^Ge-^'^'-"^) (4.4) 
and its eigenvalues are Ai^2 = e^*^ with cos A = Re A. 
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Rewriting the initial condition (^q) terms of the eigenvectors of O and substi- 
tuting in formula ( [4 .31) , one gets the number of produced fermions 

Nn=\(3n? = — ^ sin^ (nA) = sin^ ink) (4.5) 

' l-(ReA)' ^ ^ sin^A ^ ^ ^ ^ 

after the complete n-th oscillation.^^ 

We notice the presence of the envelope function 

1 _ |4|2 

E ^ , = ^ , (4.6) 

l-(ReA)^ l-(ReA)^ ^ 

which modulates the oscillating function sin^(nA). 

With increasing n this last function oscillates very rapidly and can be at all 
effects averaged to 1/2. One is thus left with the envelope function which shows the 
presence of resonance bands. The resonance bands occur where A is real and — > 1. 
It is easy to understand that their width exponentially decreases with increasing 
momenta k. To see this, let us consider the behavior of E at high momenta. In this 
regime, the function A is given by 

A^(l-e--^')e^^-, p=^^, (4.7) 



where the phase 0a can be read from eq. ([4.4]) . Near to the points where cos0a = 1 
the envelope function behaves like^^ 

E ~ - ^- - ~ = (4 8) 

1 - (1 - e-^p'f cos2 (Pa e-p' (1 - cos 0a) + 1 ' 

The width of the band can be defined as the distance between the two successive 
points at which E = 1/2. From the last expression it follows that the difference 
between the phases of A in these two points is given by Aa = 2^f2e-^'P\ Since the 
most rapidly varying term which contributes to the phase of A is (p^/2) log (p^/2), 
the width of the band can be thus estimated to be 

- -r^^ e--^ . (4.9) 
plog(p72) 

We show in figure ^ the envelope of the produced fermions in a static Universe 
for the parameters q = 10^ and mx = 100 m^. The peaks occur where A is real 
and it is confirmed that their width decreases very rapidly at increasing momenta. 
Due to the fact that the last peaks plotted are indeed very sharp, the resolution of 
figure ^ does not allow to see their top at nfc = 1. 

-'^^In ref. it is said that it is possible to extract the average over one period of the occupation 
number from the knowledge of the solution of eq. ( ^.11 ) with initial conditions u+ (t^q) — 1 1 u'^ (j]o) — 



0. The matching of this procedure with our formulae (B.2) and (4.4) gives the average time evolution 
(in our notation) N{t) = jfip sin^Kn - A)t/T]/ sin^[7r — A], where T is the period of one oscillation. 
Clearly, this result coincides with our eq. (|4.5| ) at any given complete oscillation. 
^^A completely analogous behavior occurs where cos = — 1. 
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Figure 3: Envelope of the spectrum of the produced fermions in a static Universe. The 
physical parameters are q = 10^ and mx = 100 m<^. 

5. Expansion taken into account 

As stated in the previous section, the resonance bands disappear when the expansion 
of the Universe is taken into account. In this section we will show how the occupation 
number varies when a non-vanishing Hubble parameter is considered. As we have 
seen in section ^ eqs. (|3.14|) and ( |3.15| ) give a very good agreement with the numer- 
ical results. On the other hand, the presence of phases in eq. ( p.l5| ) makes the exact 
analytical treatment of the occupation number impossible. Now, the same observa- 
tion made in the bosonic case turns very useful also to us: the phases in eq. (|3.15| ), 
when the expansion of the Universe is taken into account, are not correlated among 
themselves. As a result, the final spectra present several high frequency oscillations 
about some average function. The positions of the peaks of these oscillations depend 
on the details of the phases. However, the "mean" function can be understood in 
a surprisingly easy way. Our problem can be treated as one customary does when 
dealing with the "random walk" . Imagine one has to calculate the quantity 



S=\Ai + A2 + A3 + --- + An\ 



(5.1) 



where the Ai are complex numbers with random phases. The "random walk recipe" 
indicates that the best estimation of the above quantity is achieved by summing the 
squares of the terms Ai, since the mixed products average to zero. The chaoticity of 
the final spectra for suggests that this may also be true in our case, and this is 
confirmed by comparison with the numerical results. 
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With this method, eqs. (|3.14 ) and ( p.l5 ) turn into the much simpler relations 



_ /|a„„ 

where we remember 

\Fn\' = I - e~-P" , \Hn\^ = e--P", \Fn\^ + \Hn\^ = 1 . (5.3) 

Assuming no fermions in the initial state, and applying n times this formula, it 
is easy to see that the occupation number after the n-th production is given by 

1 1 " 

Nn{k) = -^-l[(l-2e-^^). (5.4) 

1=1 

A similar result holds for preheating of bosons, cfr. [Q where the idea of averaging 
on almost random phases is first introduced. In the bosonic case, one can exploit the 
fact that, due to the high efficiency of the production, the occupation number after 
the (n + l)-th creation is (almost) proportional to the occupation number after the 
n-th one: 

Nn+i{k) ~ (l + 2e-''")iV„(A;), (5.5) 

where the quantity Kn is analogous to our parameter p„. This simplification is not 
possible in our case, since the Pauli principle forbids Nn to be sufficiently high. 
However our final result, eq. ( |5.4| ), is also cast in a very simple and immediate form. 

The validity of eq. (|5.4] ) is confirmed by our numerical investigations, as we show 
here in one particular case. In figure ^ we compare the behavior of the "mean" 
function for the spectra with respect to the numerical one. We choose the physical 
parameters to be g = 10'^, = 4m0, and we look at the results after the 7-th 
production (this corresponds to the choices made in figure We see that the "mean" 
function interpolates very well between the maxima and the minima of the numerical 
spectrum, and that it can indeed be considered as a very good approximation of the 
actual result. 

This is confirmed by figure where we plot (for the same values adopted in 
figure ^ the quantity n^k'^ rather than the occupation number alone. This quantity 
is of more physical relevance when one is interested in the total energy transferred 
to the fermions, since the total number density of produced fermions is (apart from 
the dilution due to the expansion of the Universe, that will be considered only in the 
final result) 

Nx = ^jdkk^nk. (5.6) 

As shown in the plot, the result for Nx in the numerical and in the approximated 
case are in very good agreement. 
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After checking the vahdity of the approximation given by the "mean" function, 
we adopt it to understand how the production scales when different values of the 
parameters q and mx are considered. 

Equation ( |5.4| ) allows us to give an analytical estimate of the total amount of 
energy stored in the fermions after the ra-th production, and in particular after that 
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the whole process of non-perturbative production has been completed. Notice that 
all the dependence on the physical parameters is in the coefficients 



(7rV2p. 



(5.7) 



where k is the comoving momentum, and the only part we have to determine are the 
numbers 

.J.M„,„..,^. (5.8, 

Here, and in what follows, we express the dimensionful quantities in units of m^, the 
inflaton mass, apart from the inflaton field that is given in units of Mp. 

We also introduce the ratio R = 2q^/'^ /nix- This quantity is the most relevant, 
since it determines the zeros of the total fermionic mass and thus the values of the 
ZiS. Indeed, from eq. ( |2.2|) we see that the total mass vanishes for 0,, = —(po/R. It 
is convenient to study the production in terms of the two independent parameters q 
and R (rather than q and mx) since, at fixed R, all the spectra are the same provided 
we rescale k oc q^^^ (cfr. eq. (pl8|)). 

One can now proceed in two different ways, and we devote the next two sub- 
sections to each of them. First, one can evolve the equations for the inflaton field 
alone and find numerically the values Zi for given q and R. Inserting these values 
in eq. ( [5.4| ) one can get final values for the production which, as we have reported, 
are very close to the numerical ones. This method allows to get results which aver- 
age the actual ones, and it has the advantage of being much more rapid than a full 
numerical evolution. 

Alternatively, one can proceed with a full analytical study in order to understand 
the results given by the ffist semi-analytical method. 



5.1 Semi- analytical results 

In this subsection we evaluate the analytical formula (|5.4|), taking the coefficients Zi 



from a numerical evolution of the inflaton field. As we have said, this method gives 
results which are very close to the numerical ones, due to the fact that the "mean" 
function interpolates very well the actual spectra of the produced fermions. 

The ffist thing worth noticing is that, for each choice of q and R, the maximum Zi 
occurs at about half of the whole process of non-perturbative creation. It thus 
follows that fermions of maximal comoving momenta will be mainly produced at 
the half of the process. Our semi-analytical evaluations support this idea: the total 
energy stored in the created fermions increases slowly during the second part of the 
preheating. To see this, we show in figure ^ the numerical results for the quantity 
Nx = 2/ (ifefc^rifc/vr^ as function of the number of production. We fix the physical 
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Figure 6: Growth of Nx with the number of productions for q = 10^ and R = 30. 



parameters to be g = 10^ and R = 30, that is mx = 67 m^. With this choice, 
the total mass vanishes 10 times, and figure ^ shows the resuhs after each step of 
this production. We see indeed that the final productions are less efficient than the 
previous ones. 

From eq. (|5.4| ) it is also possible to show the evolution of the spectra with the 
number of productions. We do this in figure |^. We choose the same parameters as 
in figure ^, and we show the results after each complete oscillation of the inflaton 
field (that is after each two productions). We observe that the production rapidly 
approaches a step function in the momentum space, i.e. there exists a maximum 
momentum below which Pauli blocking is saturated (notice that the value 1/2 follows 
from the average understood in the "mean" function), and above which Uk — 0. The 
fact that the last productions do not contribute much to the total energy is also 
confirmed. 

We now turn our attention to the total energy transferred to fermions after the 
whole preheating process is completed. We fix the parameter q to the value 10^ 
and we investigate how the total integral Nx changes with different values for the 
parameter R}'^ The results are shown in figure ^ for R ranging from 5 to 10000. For 
the last value the total fermionic mass changes sign more than 3700 times and a full 
numerical evaluation would appear very problematic. This can be done in our case, 
thanks to the analytical expression ( p.4|) found, and our results extend the validity 
region of the previous full numerical study M. 



^^The behavior at small k is inessential, since this region does not significantly contribute to the 
total energy. 
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As reported, the scaling of the final result with q at fixed R is simply understood from eq. (5.S) 
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Figure 7: Evolution of the spectrum of produced fermions with the number of productions 
for q = 10^ and R = 30. 
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Figure 8: Comparison between the analytical and the semi-analytical results for Nx, at 
fixed q = l{f. 



In figure §, the results of our semi-analytical method are also compared to the 
full analytical ones of the next subsection. This comparison will be discussed below. 

From the scaling of Nx with R just reported, it is now easy to estimate the total 
energy transferred to fermions for generic values of q and R. We are interested in 
comparing our results to the numerical ones of the work [j^ . To do so, we consider the 
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Figure 9: Total final energy density produced (normalized to the inflaton one) for different 
values of q and mx ■ 



ratio between the energy density given to fermions and the one in the inflaton field^^ 

(5.9) 



We present our results in figure ^. We report them in terms of q and mx, in 
order to have an immediate comparison to the analogous figure 1 of ref. 0]. 

In figure ^ for any fixed q, the greatest plotted value for mx corresponds to 
the choice R = 5. We are not interested in extending this limit since we know 
that for greater mx (actually for values greater than the bound mx ~ v^/2) the 
production suddenly stops. The smallest value plotted for mx (at any fixed value q) 
corresponds instead to i? = 10000, that is to considering more than 3700 productions 
in the numerical evaluation of eq. ( ^.41) . 

Our final values are in good agreement with the ones of figure 1 of ref. in the 
regime of validity of the latter. The numerical results reported in that figure exhibit 
small fluctuations about an average function pxijnx)- Our results give this average 
function. This was expected, since the expression that we integrated, eq. ( |5.4| ), 
interpolates between the maxima and the minima of the numerical spectra. 



"'^^This ratio should be calculated at a time tend at the end of preheating, when the total fermionic 
mass stops vanishing. Thus in the denominator of eq. ( [3.9| ) the comoving momentum k should be 
replaced by the physical one p — k/a, with a scale factor of the Universe at tend- Analogously, the 
value c^o of the beginning of reheating should be replaced by the one at tend- However, both the 
replacements cancel out in the ratio, since both and redshift as energy densities of matter. 
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The numerical results of ref. have a smaller range of validity than the ones here 
presented. This occurs because the numerical evolution of that work is limited to the 
first 20 oscillations of the infiaton field, and so fermionic production has not come to 
its end for small values of rrix- We confirm that at high values of mx (actually at 
small values of R for any given q) the production depends very weakly on mx- In 
addition, our results show a decrease of the energy transferred to fermions for smaller 
values of mx- This behavior will be explained in details in the next subsection. 

5.2 Analytical results 

We want now to show that all the results presented in the previous section can be 
also achieved with a full analytical study of eq. ( |5.4| ). 

First of all, we have to estimate the quantities Zi given in eq. ( ^.8|) . To do 
this, it is more convenient to work in terms of the physical time t: after the first 
few oscillations, the infiaton evolution is very well approximated by the expression 
(remember t is expressed in units of m'^^, while in units of Mp) 

, / N 1 cos (t) ,^ ^ 

0(t)^-=— ^. (5.10) 



The scale factor of the Universe follows the "matter-domination" law, and it is well 
approximated by a(t) = t}!"^ . 

The values t*j are determined by the condition of vanishing of the total mass of 
the fermions, that is, by making use of eq. ( p.lQ ), 

cos(t,i) = (5.11) 



where we remind R = (2^) /mx- The parameter A = (yovr^o)"^ is of order one 
and will not play any special role in what follows. Notice that the last production 
occurs at t ~ RA. 

Hence, keeping only the dominant contribution to the derivative of (p with respect 
to the physical time, we get the expression 

.i.M«v3^v3(|^y'yi_(l^y, (5,2, 

where we can assume t^i ~ in. 

Equation (|5.12|) exhibits a very good agreement with the numerical evaluation 
of the same quantity. It also shows that the maximal value for Zi is reached at = 
AR/2, that is, at half of the whole process of non-perturbative creation. This was 
anticipated in the previous section, where we showed that the most of the fermionic 
production occurs in the first part of preheating. 
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Starting from eq. ( ^.12| ) we can also calculate the number density of produced 
particles 



Nx{q,R) 



dkk^N{k) 



(5.13) 



where N{k) is obtained from eq. ( |5.4|) with n = nmax = {RA)/7t. 

For R big enough, the product in eq. (|5.4|) can be written as the exponential of 
an integral. Thus, we obtain 

'9 /1 4/3 \ 3/2 

-j q'/'R'/'x 



N- 



mx 



X 



47r3 \^ TT 



(5.14) 









exp 


/ dy log 


1 - 2 exp 1^ 









"^1/3^/13^ 



with the substitution t = k- ^/V/{2J^lhp^^^ . 

The integral in dy which appears in eq. ( ^.14| ) cannot be calculated analytically. 
Anyhow, we can approximate it by 







/ dy log 


1 - 2 exp 1^ 


/o 





^1/3^/13^2 



9 {t) lo^ 



2e 



-1.5i^ 



(5.15) 



where g{t) is a function of order one that, for our purposes, can be approximated by 
a constant c in the range 0.5 < c < 1. 

Hence, the integrand within curly brackets in eq. ( p.l4| ) rewrites 



1 



1 - 2e 



(5.16) 



In the large- i? limit this function approximates a step function, which evaluates to 
one for 



TT log 2 

2RAc 



< t < 



(5.17) 



and to zero for the remaining values of t. 

Since the quantity ( ^.16| ) is proportional to the occupation number rik, our an- 
alytical calculation confirms the usual assumption that, after few oscillations, the 
fermionic production saturates the Fermi sphere up to a given maximum momentum 
^max- This is also shown in the previous figure |^. However, our derivation gives a 
different scaling for fcmax with respect to the previous hterature |jl2|, 0. Indeed, from 
eq. (|5.17| ) it follows (apart from proportionality factors) 



,1/3 



m 



1/6 
X 



'log 



,1/2 



mx 



(5.18) 



which is however quite close to the result given in [0 . We will comment more about 
the origin of the scaling (|5.18| ) in the conclusions. 
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From eqs. ( |5.14 ) and ( 5.17 ), we obtain the final expression for the number density 
of the fermions created during the whole process 



Nx{q,R) 



3^ V l.Svr ) 



3/2 



/2 



log 



1/2 



AAc q 
n log 2 nix 



-\ 3/2 



(5.19) 



We can now go back to figure ||, where this last equation (called "analytical" in 
the figure) is compared to the result with the semi-analytical method of the previous 
subsection. In plotting eq. (|5.19|) we chose c = 0.78 for the numerical factor involved. 
As we can see, the final results achieved with the two methods are in very good 
agreement with each other, thus confirming the validity of the formula (|5.19|) .^^ 

Rewriting eq. (|5.19|) in terms of q and mx we can draw some conclusions. First, 
apart from a logarithmic correction, the scaling of the total energy 



1/2 



px oc nixNx oc qnix 



log 



iAc q 



1/2 



7rlog2 mx 



3/2 



(5.20) 



is linear in q, as generally expected |^ . The dependence of px on mx requires some 
more care: the threshold value for mx is given by the condition ~ 4, that is,^'' 



(mx) 



th 



Vq 

2 



(5.21) 



For values of mx not much smaller than {mx)th, the total energy depends very 
weakly on mx, this result being in agreement with the numerical evaluation given 
in 1^. On the other hand, for values of mx much smaller than (mx)th, the factor 



1 /2 

m X starts to dominate, and we expect it to determine the scaling of the total energy 
when mx — > 0. 



6. Backreaction 

The results presented so far have been achieved neglecting the backreaction of the 
produced fermions on the evolution of the infiaton field and on the scale factor. 
This is a common approximation, since a more complete treatment (especially an 
analytical one) of the whole phenomenon is a very difficult task. However, the effects 
of the backreaction can be understood to a good degree of accuracy in the Hartree 
approximation. 

For what concerns preheating of fermions, the Hartree approximation consists in 
taking into account the term 

9{XX) (6.1) 

^^The small discrepancy between the two curves can be attributed to the fact that c is not exactly 
constant. 

^''The number 4 comes from the fact that the value of the infiaton at its first minimum is (/) ~ 
-0.07Mp, while at beginning 0o = 0.28Mp. 
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into the evolution equations for the inflaton and the scale factor. The equation for 
the field thus rewrites (in physical time) 



+ 3ff0 + m20 + ^(XX) = 0. (6.2) 

The study of this effect has been performed numerically in ref. 0, where it is 
shown that backreaction starts to be important for q > qf, 10^ — 10^°. Figure 7 of 
that work shows how the evolution of the inflaton field (j) is modified when backreac- 
tion is considered and q is sufficiently high. First, one observes that the amplitude 
of the oscillations of the infiaton is very damped already after the first production. 
This effect is the most obvious one, since the term (|6.1| ) takes into account the decay 
of the inflaton into fermion-antifermion pairs, while in its absence the equation for 
considers only the damp due to the expansion of the Universe. The second feature 
that emerges from the evolution performed in ref. is that at the beginning the 
fleld does not oscillate about the minimum of the potential V = m|0^/2, but about 
the point 0* where the total fermionic mass vanishes. Moreover, the frequency of 
these oscillations is higher than m^. 

These last two effects are due to the change in the effective potential for 
induced by the term (|6.1|), and disappear when the quantity {XX) is decreased by 



the expansion of the Universe. Their net effect is to render the whole mechanism 
of preheating more efficient, since the rise in the frequency of the oscillations of 
increases the number of productions of fermions. In ref. [0 it is indeed shown that for 
q = 10^ the total production is about 5% larger than the one without backreaction, 
while for q = 10^° the increase is about 50%. 

We now briefiy study the evolution of the infiaton under eq. (|6.2|) by means 
of the analytical results presented above. We show that even a very approximate 
analysis confirms the numerical result that indicates in g ~ 10^ — 10^° the threshold 
above which backreaction should be considered. 

To begin on, the term ( |6.1| ) needs to be normal ordered. Doing so, one gets |^ 

(XX) = fd^]^(l + — - . (6.3) 



{2TTaf J \ ' to 



In terms of the Bogolyubov coefficients, this quantity evaluates to 

+ Re(^«/?^e2^/-'^^^ 



{XX) = —-3 / d'k 
[ZTTaj 



2 am 

to 



(6.4) 



We see that {XX) vanishes for f3 = 0. This is obvious, since backreaction starts 
only after fermions are produced. Some approximations can render eq. (|6.4| ) more 
manageable. First, we notice that the oscillating term in the exponential averages 
to very small values the integral of the second term in square brackets. Second, 
we see that k <^ |am| where jSk is significantly different from zero. From both 
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these considerations, the integrand in eq. ( |6.4|) can be approximated (up to the sign 
of m) by the occupation number |/3p, so that the whole effect is (approximatively) 
proportional to the number of produced fermions. 

The numerical results of ref. show that, for the values of q for which the 
backreaction is to be considered, its effects can be seen already in the first oscillation 
of the inflaton field. Since we are only interested in estimating the order of magnitude 
of Qb, we thus concentrate on the first oscillation of 0, neglecting the expansion of 
the Universe in this short interval. 

With all these approximations, eq. (|6.2|) rewrites 



2/" + y + 10-iV/^A = 0, (6.5) 
\m\ 

where we have rescaled y = (p/cpQ and we remind that the time is given in units of m'^^. 

This equation is very similar to the one obtained in the bosonic case The last 
term changes sign each time m = and, when sufficiently high, forces the inflaton 
field to oscillate about the point 0,, at which the total fermionic mass vanishes. It 
is also responsible for the increase of the frequency of the oscillations. To see this, 
we assume that this term dominates over the second one of eq. ( |6.5|) and we solve 
it right after the first fermionic production at the time In the absence of the 
second term, eq. ( |6.5| ) is obviously solved by a segment of parabola until the time 

— + 2|?/'(t*i)|/(g^/^10~"'^^) at which m vanishes again and the last term of 
eq. ( |6.5| ) changes sign. As long as the second term of eq. ( p.5|) can be neglected, the 
inflaton evolution proceeds along segments of parabola among successive zeros of the 
mass m. The "time" duration of these segments is expected to be of the same order 
of the first one, since the successive fermionic productions balance the decrease of 
{XX) due to the expansion of the Universe. 

The period of these oscillations can thus be roughly estimated to be T ~ 
2(t*2 ~ We see that, for q > 10^, this period is smaller than the one that 

the inflaton oscillations would have neglecting backreaction. Since this increase of 
the frequency is the main responsible for the higher fermionic production, the result 
Qb ~ 10^ can be considered our estimate for the value of q above which backreaction 
should be taken into account. 

This result, although obtained with several approximations, is in agreement with 
the numerical one of ref. [0 . 



7. Conclusions 

Preheating of fermions and bosons present several analogies. In both cases, especially 
when one is interested in the production of very massive particles, the creation occurs 

^^However, after the first part of the process, the production looses its efficiency and the expansion 
of the Universe dominates. As we have said, the term {XX) can then be neglected and the inflaton 
starts oscillating about the minimum of the tree level potential. 
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for very short intervals of time, during which the frequency of the particles varies non 
adiabatically and their occupation number cannot be defined. This strong similarity 
can be exploited to extend to the fermionic case the formalism developed in ref. for 
the analytical study of preheating of bosons. This was first done in ref. [^|, where 
the result for the particles created in the first production is reported. When one 
wants to consider the successive productions, a more detailed analysis is necessary. 
This has been the object of the first part of the present work. The formulae that we 
obtained, valid for an arbitrary number of productions, give a very good agreement 
with the numerical results, as some examples provided manifest. 

Another strong similarity between the two phenomena is that in both cases the 
production would occur only through resonance bands (in momentum space) was 
not for the expansion of the Universe. For what concerns fermions, this feature was 
first studied in ref. |Tl], |T2|, where production in a static Universe is considered. As 
it is well known, the expansion of the Universe removes the resonance bands and, 
as a result, this allows the fermionic production to saturate the whole Fermi sphere 
up to a maximal momentum fcmax- Our analytical results confirm both the presence 
of the resonance bands in the static case (where our formulae considerably simplify) 
and their disappearance when the expansion is considered. 

Rather than the detailed form of the spectrum, a quantity which may be of more 
physical interest is the total energy of the produced particles. The analogy with 
the bosonic case turns useful also in this regard. In the final spectra one notices 
the presence of some maxima and minima, whose exact position is determined by 
the detailed knowledge of the phases of the Bogolyubov coefficients (adopted in the 
analytical derivation) after each production. However, due to the expansion of the 
Universe, these phases are effectively uncorrelated among themselves. This gives the 
production a stochastic character which is responsible for the disappearance of the 
resonance bands. Averaging over these phases, it is possible to get an analytical 
function for the mean occupation number. For what concerns preheating of bosons, 
the results are simplified by the fact that, in the region of very efficient production, 
the mean occupation number after a given production is proportional to the mean 
occupation number of the previous one. This assumption is not possible in the 
fermionic case, since the Pauli principle prevents the occupation number to exceed 
1. However, it is remarkable that also for fermions the final result can be cast in a 
very simple and readable form. 

The analytical formula for the mean occupation number confirms the saturation 
of the Fermi sphere (times the factor 1/2 which comes from the average) up to a 
momentum fcmax- In order to get this quantity, it is sufficient to know the derivative of 
the infiaton field and the value of the scale factor at the points where the productions 
occur. These values can be calculated numerically evolving the equations for the 
infiaton field alone (at least when backreaction is neglected), or analytically. The 
former possibility is much more rapid than a full numerical computation, since one 
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does not have to consider the equations for the fermions. Moreover, its results agree 
very well with the full numerical ones 0. Thanks to the increase in the rapidity, 
this evaluation (that is called "semi-analytical" in the present work) allows to get 
results with a more extended range of validity with respect to the previous full 
numerical ones. 

The results provided in this way can also be achieved with a full analytical study. 
From both these methods, we deduce that the ratio px / Pd> scales as 



PX 1/2 



log 

mx 



3/2 

(7.1) 



up to (mx)th ~ 

For nix not too smaller than {mx)th ^^is scaling is in good agreement with 
the one px/ Ptf, ^ q oi the numerical work |]^. However, we see that the density of 
produced fermions decreases at smaller values of mx- 

We would like to conclude our work showing that the scaling ( |7. 1| ) can be also 
achieved from very immediate considerations. As we said, preheating of fermions 
saturates a Fermi sphere in momentum space up to a maximal momentum fcmax- 
From the analytical formula ( |5.2|) , we notice that at high momenta k the occupation 
number is well approximated by 

n 

iV„(A;)^5^e-^"/^^ (7.2) 

i=l 

where we remember Zi oc q^^^a}^'^ {ri^i)\(f)' [ri^i]]^/'^ . 

In this last equation, we replace all the parameters Zi with a mean value z, so 
that Nn ~ nexp( — k'^/z^). The scaling of z with the physical parameters q and mx 
follows from the scaling of all the Zi. The maximal momentum k^^^ is thus expected 
to scale as the quantity 2;j(logn)^/^. Considering now the evolution of the infiaton 
field in physical time t, we notice that both the number n of productions and the 
times t*j at which they occur are proportional to the parameter R = g^/^/(2mx). 
Moreover, we see that the Zj^s scale as 



Zi oc g^'^'^a^, 



dt 



''^'ocg^/^tf47 0cgV4^V6. (7.3) 



3 



max J 



We thus get fcmax OC g^''^i?^/^[log-R]^/^, from which (remember px — mxk 
the scaling ( [7.1| ) simply follows. 
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